Percolation-dependent Reaction Rates in the Etching of Disordered Solids 



K. M. Kolwankar 1 - 2 , M. Plapp 2 and B. Sapoval 1,2 
1 Centre de Mathematiques et de leurs Applications, Ecole Normale Superieure - CNRS, 94-14-0 Cachan, France 
2 Laboratoire de Physique de la Matiere Condensee, Ecole Polytechnique - CNRS, 91128 Palaiseau, France 



(N 
O 
O 
(N 

> 
O 



> 

OS 
(N 



(N 
O 



c 

o 
o 



A prototype statistical model for the etching of a random 
solid is investigated in order to assess the influence of disorder 
and temperature on the dissolution kinetics. At low temper- 
ature, the kinetics is dominated by percolation phenomena, 
and the percolation threshold determines the global reaction 
time. At high temperature, the fluctuations of the reaction 
rate are Gaussian, whereas at low temperature they exhibit a 
power law tail due to chemical avalanches. This is an exam- 
ple where microscopic disorder directly induces non-classical 
chemical kinetics. 
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Chemical etching of disordered solids is a process of 
great practical importance pi. In some cases, such as in 
corrosion, it can be a cause of undesirable deterioration. 
In other cases it is useful, for example when one wishes to 
roughen a smooth surface for adhesion or to obtain some 
specific optical or esthetic properties. Such processes are 
also stimulating from a theoretical point of view since it 
has been observed in experiments that they can 

give rise to complex patterns and dynamics. 

To understand the universal features of random etch- 
ing processes, the study of simplified toy models is useful 
since they allow to grasp some essential results using the 
concepts of percolation theory ||. The randomness of 
the solid may be caused either by the topology of the 
chemical bonds, such as in amorphous media, or by fluc- 
tuations of the composition, such as in a solid solution. 
It can also be the randomness of a passivation layer that 
may appear spontaneously in the corrosion of polycrys- 
talline solids. 

In this letter, we focus on the time evolution of the 
etching process and the influence of temperature. In 
the previous studies of random etching, random site re- 
sistances were assigned, and "strong" sites were never 
etched §,@,|. This has given a theoretical interpreta- 
tion to the spontaneous stabilization of fractal surfaces 
in a model where the etchant is consumed in the course of 
the reaction |^,0j8) . The fluctuations observed in || were 
interpreted in terms of a self-organized fractal growth 
linked to percolation theory j7[] . 

However, if the disorder is at some microscopic level, 
it is more realistic to assign random reaction rates to 
the dispersed random elements of the solid. Then, even 
a site with a very small dissolution rate will eventually 
be etched. Furthermore, since real reaction rates are 
temperature-dependent, the influence of this parameter 



can be studied. In the following, we find a transition from 
smooth and continuous growth at high temperatures to 
a highly intermittent regime at low temperatures that 
exhibits chemical avalanches. 

The 2d square lattice model of the solid consists of 
elements i that have a random site-dependent binding 
energy Ei. The dissolution rate Ri of a molecule i in 
contact with the etchant is supposed to follow an Arrhc- 
nius law that corresponds to an activated process: 



Ri = C lu cx_-p(-Ei/k B T) 



(1) 



where C is the concentration of the etchant, ks is Boltz- 
mann's constant, and T is the temperature. For sim- 
plicity, the attempt frequency loq is assumed to be the 
same for all sites. The energies Ei are supposed to be 
distributed over a finite range. The simplest choice is 
Ei = Eq +piAE where the p^s are distributed uniformly 
between and 1. Then, Eq. [l] can be rewritten as 



Ri = C lo cxp(-EpPi) , 



(2) 



where uo — exp(— Eq/HbT) and Ep = AE/ksT. Note 
that this simplified model neglects many factors that may 
affect real dissolution rates, such as for example steric 
blocking. There are essentially two limits depending on 
the value of Ep. With ksT of order 25 meV at room 
temperature, weak disorder ("high temperature") typi- 
cally corresponds to Ep of order unity, and strong chem- 
ical disorder ("low temperature") corresponds to Ep of 
order 10 to 100. 

We consider an irreversible reaction with no redepo- 
sition, and we suppose that one molecule of etchant is 
consumed per dissolved site. It follows that for a finite 
volume of etchant with initially N etchant molecules and 
concentration Co, the concentration after the dissolution 
of j sites is Cj = C*o(l — j/N). It therefore decreases 
with time, and the reaction stops when all the etchant 
molecules are exhausted. Furthermore, we suppose that 
the mixing in the solvent (by diffusion of the etchant 
or hydrodynamic flow) is fast compared to the chemi- 
cal events, such that the concentration in the solution 
remains uniform. 

The time evolution of etching is studied through an 
event-oriented Monte Carlo algorithm in physical time. 
A site i dissolves within a time interval At with probabil- 
ity RiAt if At is small enough (RiAt <C 1). Since these 
probabilities are independent for different sites, the prob- 
ability that any of the surface sites is etched in the time 
interval At is X)ig5 ^i^i where the sum runs over all the 
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FIG. 1. Snapshot pictures of the etching front at the end 
of the reaction for TV = 10 5 , L — 10 3 and (a) Ep = 1, (b) 
Ep = 10 and (c) Ep = 100. The unetched solid is black 
and the solution is white. Note that the surface of the white 
regions is the same and equal to N, the initial number of 
etchant molecules in the solution. 



surface particles (active sites). The time interval before 
a site is etched with probability one is therefore 



At 
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1 



J2ies CujJ q h(p)exp(—Epp)dp 



(3) 



where h(p) gives the number of active sites with resis- 
tance between p and p + dp. Within that time interval, 
a given site i is etched with probability 



Pi = 



Rj 



exp(—Eppi) 



J2ies Pi Jq h(p) exp(-Epp)dp 



(4) 



In each Monte Carlo step, the algorithm chooses to etch 
a surface site i with probability Pj and increases the time 
t by At. As a consequence, the physical time is not pro- 
portional to the number of Monte Carlo steps, since At 
depends on the dissolution rates of all the surface sites. 
However, by construction, the average time needed to 
dissolve a given site with reaction rate Ri is always the 
same, independently of the system size or the surface 
configuration. The probability distribution of Eq. ^ has 
already been used in Ref. || to study the influence of fi- 
nite temperature on invasion percolation, but this study 
was restricted to the initial stages of the process. 

The initial lattice of width L has a flat interface. The 
ratio N/L is the etching depth if uniform dissolution oc- 
curs. We study here the case N/L ^> 1 such that deep 
pores can be formed. In this situation, the essential pa- 
rameter that governs the pattern formation is the reduced 
energy Ep which controls the relative dissolution speed 
of "hard" and "weak" sites. For small Ep, the dissolu- 
tion rate is nearly the same for all sites whatever their 
energy; in the limit Ep — > 0, this model maps onto the 
Eden model M]. In this regime, etching leads to the for- 
mation of a microscopically rough, but macroscopically 
smooth interface that has a finite width W and prop- 
agates with a roughly constant speed (Fig. [l]a). When 
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FIG. 2. Activation plot of the reaction half-time T\/2- 

Ep increases, sites with small Pi become relatively much 
faster than sites with large p,-. The dynamics leads to 
a irregular or porous interface (Figs. and |l]c), with 
disconnected islands surrounded by sites which may sur- 
vive for a long time. These islands continue to be slowly 
etched and their size therefore decreases with time. In 
the limit Ep — > oo, the weakest site is chosen with prob- 
ability 1 and the model becomes identical to classical 
invasion percolation fll| . 

For this reason, it is not surprising that the percolation 
properties of the underlying lattice enter the dynamics in 
the large Ep limit. Our main focus here is the chemical 
kinetics. The global reaction rate strongly depends on 
temperature. To characterize its behavior, we show in 
Fig. H the reaction half-time T1/2, that is, the time it 
takes to consume half of the N etchant molecules, as a 
function of Ep. For large enough Ep, this plot gives an 
cffcctice activation energy equal to (0.59 ± 0.01), close to 
the square lattice percolation threshold p c . The time for 
the reaction to go from N/2 to N/4 etchant molecules 
follows the same behavior. Simulations for the 3D cubic 
lattice give a similar result with a slope 0.30 ±0.01, close 
to p c for this 3D lattice. This indicates that, for strong 
disorder (or low temperature) the global reaction has an 
activation energy equal to Eq + p c AE (the Eq contribu- 
tion comes from the definition of u>). To our knowledge, 
this is the first direct appearance of a percolation thresh- 
old in the time domain. 

Figure |^ shows the histograms h(p) of the site resis- 
tances pi on the solid-liquid interface (including the is- 
lands) after 50% of the etchant is consumed for different 
temperatures. They count the number of surface sites 
with p < Pi < p + e as a function of p (here e = 0.01). 
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FIG. 3. Histogram of piS on the interface after 50% of the 
etchant is consumed (averaged over 10 realisations). Here, 
N = 10 6 , L = 10 3 . 

For small Ep = 1 (high temperature), h(p) increases with 
p, but there is always a large number of weak sites on 
the surface. In contrast, for Ep 3> 1 (low temperature), 
there are very few weak sites, and h(p) starts to increase 
sharply at a value close to the percolation threshold. 

These data help to explain the appearance of the per- 
colation threshold in the reaction rate. The total time 
needed to etch n molecules is simply the sum of the time- 
steps for each individual etching event, t — Ylj=i^>j- 
At high temperature, all the dissolution rates are of sim- 
ilar magnitude. Therefore, the time-step depends only 
weakly on the surface configuration and the reaction pro- 
ceeds with a well-defined average rate. In contrast, at low 
temperature, the contributions of successive configura- 
tions are very different. For a typical configuration that 
contains a few weak sites, the sum in the denominator of 
Eq. (||) is dominated by the large rates associated with 
these sites, the corresponding time-step is small, and the 
weakest site is etched with a probability close to one. 
However, from time to time, the growth front coincides 
with the perimeter of a percolation cluster. In this situa- 
tion, the weakest site has a strength close to p c , and the 
corresponding time-step scales as At ~ exp(_E^p c ). As 
a consequence, the dynamics consists of series of rapid 
etching events ( "chemical avalanches" ) separated by long 
time steps when the front has reached the surface of a 
percolation cluster. The appearance of the percolation 
threshold in the reaction time indicates that the latter 
are the rate-limiting steps. 

The time evolution of the histograms can be qualita- 
tively understood on a coarse-grained time scale, that is, 
a time scale much larger than individual time-steps. On 
this scale, one can define a global reaction rate R, that is, 
the total number of sites etched per unit time. Note that 
R must be proportional to Civ. One can then write a 
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FIG. 4. Histogram of piS on the final interface for Ep = 10, 
TV = 10 6 and L = 10 3 . 

continuous evolution equation for the number of surface 
sites of strength p: 

^ = Re - h(p)Cive X p(-Epp). (5) 

The first term counts the arrival of freshly uncovered sites 
(with uniform distribution) on the surface, whereas the 
second term is the product of the etching rate of sites 
with strength p times their number. One clearly sees 
that h(p) = Reexp(Epp)/(ivC) is a steady-state solution 
of this equation. This is indeed what is observed for 
high temperatures. For low temperatures, inserting R ~ 
Civexp(-E p c ) yields h(p) ~ eexp[Ep(p - p c )\ @, in 
remarkable agreement with the final histogram data for 
Ep = 10, shown in Fig. [|. This has a general consequence 
that could be called the "chemical survival of the fittest" : 
since almost all weak sites have been removed from the 
surface, the disordered solid is stronger once etched than 
it was before the attack. This has already been discussed 
in a simpler context || . 

One notes that, for even lower temperatures, the tran- 
sitory h(p) shown in Figure [| levels off for large p. This 
occurs because in our simulations the total number of 
steps is far too small to reach the steady state at low 
temperature. Indeed, for such a steady state to appear, 
a sufficient contrast must occur between sites above p c . 
For this to happen, a huge number of strong sites must be 
uncovered, and this requires an exceedingly long time. At 
zero temperature, which is equivalent to invasion perco- 
lation, the steady state is never reached. The plateau for 
the highest p values in the transitory histograms simply 
reflects the initial uniform distribution of site strengths. 

The chemical avalanches that occur at low tempera- 
tures can be characterized in two different ways. The 
first is purely geometric. When the "hard" perimeter of 
a percolation cluster is broken, the reaction has a large 
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FIG. 5. Statistical distribution of the reaction rate for dif- 
ferent temperatures: distribution of the number n of sites 
etched in a time interval t% - ti = T iot /1000. Here N = 10 6 , 
L = 10 3 . The first 20% of events have been discarded. 



experiment s p] and were explained in terms of a sand- 
pile model |l4[ |. Our findings might offer an alternative 
interpretation for these observations. 

In summary, we have introduced a minimal model for 
the etching of a disordered solid. We have studied the 
influence of temperature on chemical kinetics and etch- 
ing patterns. At high temperature (or weak disorder), 
the reaction rate exhibits gaussian fluctuations around 
a well-defined average rate, and the dissolution front is 
macroscopically smooth. At low temperature (or strong 
disorder), the reaction rate has strong fluctuations caused 
by chemical avalanches. The global reaction rate is shown 
to depend on temperature with an activation energy di- 
rectly related to the percolation threshold. Interestingly, 
the avalanches, which may be large at low temperature, 
do not influence the reaction rate itself. The effect of 
possible redeposition and that of concentration gradients 
on reaction kinetics and patterns should be examined in 
the future. 



probability to continue at the place where the last etch- 
ing event has occurred since new weak sites may be un- 
covered. The avalanche size can then be determined by 
"connected cluster events" . Starting with a given etch- 
ing event, if the next etched site is a first neighbor, it 
is added to the cluster, and so on, until the next etched 
site is not a neighbor to any of the cluster sites. For large 
values of Ep, the number of avalanches versus their size 
follows a power law with a slope around 1.53 ± 0.01 (not 
shown here), in agreement with the result for classical 
invasion percolation |jl3f . 

A more physical way to characterize avalanches is via 
the global reaction rate. It can be written as 



R 



n(t 2 - h) 



(6) 



where n(t 2 — ti) is the number of sites etched between 
times t\ and t 2 . So defined, this quantity is proportional 
to an experimental signal, for example the corrosion cur- 
rent if the etching is the result of an electrochemical pro- 
cess. To study the behavior of the reaction rate, we have 
performed simulations with a constant etchant concen- 
tration C (corresponding to an infinite volume of solu- 
tion) up to a fixed number of etched sites N . The total 
reaction time T tot was divided in equal time intervals 
< 2 — <i = T tot /1000, and the number of sites etched in 
each interval was counted (discarding the initial etching 
of the flat surface). The resulting histograms are shown 
in Fig. ||. For small E@ values, the reaction rate has 
a Gaussian distribution around a well-defined average. 
In contrast, for very large values of Ep, the histogram 
becomes close to a power law due to the occurence of 
avalanches with a broad size distribution. Therefore, the 
reaction rate has large fluctuations in both space and 
time. Recently, chemical avalanches have been found in 
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